version 10.0
log using "C:\fsu2010\formateur\replication\stata_mxl_predictprob.log", replace
#delimit ;
set mem 30m;
clear matrix;
clear mata;

set more off;

*     ***************************************************************** *;
*     ***************************************************************** *;
*       File-Name:  stata_mxl_predictprob.do                            *;
*       Date:       November 9, 2010                                    *;
*       Author:     GG/MG                                               *;
*       Purpose:    Produces Figure 1a in AJPS paper                    *;
*       Input File:     europe.dta                                      *;
*       Output File:    stata_mxl_predictprob.log                       *;
*       Data Output:    none                                            *;
*       Previous file:  europe.dta                                      *;
*       Machine:        Matt's desktop                                  *;
*     ****************************************************************  *;
*     ****************************************************************  *;

use "C:\fsu2010\formateur\replication\europe.dta", clear;

*     ****************************************************************  *;
*       Model 6: Western Europe                                         *;
*     ****************************************************************  *;

mixlogit choice largestparty_east party_seatshare_east medianparty1_east presidentparty_east previousPM_east previousPM_conflict_east previousPM_gain1_east previous_cabinet_gain1_east,  
    rand (largestparty party_seatshare medianparty1 presidentparty previousPM president_investiture president_direct previousPM_conflict previousPM_gain1 previous_cabinet_gain1) 
    ln(0) nrep(125) group(cabinetcode);
    
*     ****************************************************************  *;
*       Now produce Figure 1                                            *;
*     ****************************************************************  *;

matrix coeffs = e(b);
matrix covmat = e(V);

matrix bdraw = coeffs;
matrix bdraw[1,19] = 0;

keep if cabinetcode==1222;
set obs 10000;

generate pvda = 0;
replace pvda = 1 if partyID==7;
generate cda = 0;
replace cda = 1 if partyID==13;
generate vvd = 0;
replace vvd = 1 if partyID==22;
generate d66 = 0;
replace d66 = 1 if partyID==10;

generate ppvda = .;
generate pcda = .;
generate pvvd = .;
generate pd66 = .;

local a=1;

while `a'<=1853 {;

  matrix bdraw[1,9] = ((`a'-716)/100);
  mixlpred2 pprobs; 

  qui generate pvdat = pprobs * pvda;
  qui generate cdat = pprobs * cda; 
  qui generate vvdt = pprobs * vvd;
  qui generate d66t = pprobs * d66;

  egen pvdat2 = sum(pvdat);
  egen cdat2 = sum(cdat);
  egen vvdt2 = sum(vvdt);
  egen d66t2 = sum(d66t);

  qui replace ppvda = pvdat2 if _n==`a';
  qui replace pcda = cdat2 if _n==`a';
  qui replace pvvd = vvdt2 if _n==`a';
  qui replace pd66 = d66t2 if _n==`a';

  drop pvdat cdat vvdt d66t
       pvdat2 cdat2 vvdt2 d66t2 
       pprobs;
    
  local a = `a' + 1;
  display "." _c;

    };

drop if _n>1853;

gen MV= ((_n-716)/100); 
    
graph twoway line ppvda MV, clpattern(solid) clwidth(medium) clcolor(black) 
    || line pcda MV, clpattern(longdash_dot) clwidth(medium) clcolor(black) 
    || line pvvd MV, clpattern(longdash) clwidth(medium) clcolor(black) 
    || line pd66 MV, clpattern(shortdash) clwidth(medium) clcolor(black) 
    ||, 
    xlabel(-8 -6  -4  -2  0 2  4 6  8  10 12, nogrid) 
    ylabel (0 0.2 0.4 0.6 0.8 1)  
    legend(off)  
    yline(0 0.2 0.4 0.6 0.8 1, lcolor(white)) 
    xtitle("" , size(2.5))  scheme(s2mono) graphregion(fcolor(white)); 
    
log close;



  


